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Heavy carbon doping of MgB2 is studied by first principles electronic structure studies of two 
types, an ordered supercell (Mg(Bi_ a; C :E )2, £=0.0833) and also the coherent potential approximation 
method that incorporates effects of B-C disorder. For the ordered model, the twofold degenerate 
CT-bands that are the basis of the high temperature superconductivity are split by 60 meV (i.e. 7 
meV/% C) and the a Fermi cylinders contain 0.070 holes/cell, compared to 0.11 for MgB2. A virtual 
crystal treatment tends to overestimate the rate at which a holes are filled by substitutional carbon. 
The coherent potential approximation (CPA) calculations give the same rate of band filling as the 
supercell method. The occupied local density of states of C is almost identical to that of B in the 
upper 2 eV of the valence bands, but in the range -8 eV to -2 eV, C has a considerably larger density 
of states. The calculations indicate that the a Fermi surface cylinders pinch off at the zone center 
only above the maximum C concentration x ~ 0.10. These results indicate that Mg(Bi_ a: C :E )2 as 
well as Mgi-a, Ala;B2 is a good system in which to study the evolution of the unusual electron-phonon 
coupling character and strength as the crucial a hole states are filled. 
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I. INTRODUCTION 

Discovery of superconductivity at T c ks 40 K in 
MgB2 has enlivened not only interest in new classes 
of superconductors with high T c and novel two-gap 
behavior, but also pursuit of new materials for ap- 
plications that require high critical current densi- 
ties and high critical fields (H c2 ). In highly resis- 
tive films (due to unreacted components and/or oxy- 
gen and carbon impurities) the perpendicular and 
parallel critical fields have reached^ « 34 T, 
H^ 2 w 49 T, making MgB 2 a real possibility for a 
high field conductor. Such application require de- 
fects or grain boundaries to pin the vortices whose 
motion would otherwise lead to energy dissipation 
and joule heating. Such defects also affect the un- 
derlying electronic structure and pairing interaction, 
at least in some parts of the sample, and often de- 
crease T c . The simplest defect to understand is the 
substitutional impurity, which provides a means of 
varying the intrinsic properties in a continuous and 
controllable manner and may enhance pinning mech- 
anisms. Although readily synthesizablc as a stoichio- 
metric, rather clean compound, the MgB2 lattice re- 
sists most attempts to alloy, both on the Mg and B 
sites, at more than the very dilute level. 

The two established exceptions are Al substitu- 
tion for Mg, and C substitution for B, both of which 
lead to a rapid decrease in T c . Published reports 
on the behavior vary considerably with synthesis 
method and sample treatment, but some of the gen- 
eral behavior seems to be established. Al substi- 
tution should be the simpler one, since the very 
strongly bonded honeycomb structure B layers re- 
main intact, and simple rigid band filling of the hole 
states at first appears to be a likely possibility. Up 
to y = 0.10 in Mgi_ a Al a B2 that appears to be the 



case, with T c decreasing smoothly. Beyond the con- 
centration y = 0.10, however, there are signs of two- 
phase behavior (two c lattice parameters in diffrac- 
tion studies) arise^ Around y = 0.25 the system re- 
verts to a single phase very low and finally vanishing 
T c . A study of the trends in electron-phonon cou- 
pling strength using first principles calculations in 
the virtual crystal approximation 3 reported a sharp 
change of behavior at y— 0.25. A thorough study 
of the energetics of Al substitution revealed a strong 
tendency for superstructure formation^ at y=0.25 
(and y=0.75). Still, the observed onset of two- 
phase behavior already at ?/=0.10 is unexplained, 
and nonstoichiometry 5 and microscopic defects must 
be kept in mind. 

Unlike the Al alloying case where there has been 
at least rough consensus on changes of properties, re- 
ports of the change in T c and structure with addition 
of carbon on MgB2 have varied widely. These differ- 
ences seem in some cases to reflect real differences 
in materials due to the various methods of synthesis 
and heat treatment, and the samples prepared were 
highly mixed phase, thereby creating additional un- 
certainty. The very similar xray scattering strength 
of the B and C atoms has rendered standard xray 
diffraction ineffective in determining the C content 
of a sample. Reports early on suggested C misci- 
bilitics in Mg(Bi_ 2: C 2: )2 as smah— as x—0.02 to as 
largai as x=0.3. Recent studies of several groups 
suggest that £=0.2 can be achieved while larger con- 
centrations are questionable 8 910 although reports of 
larger C concentrations under high pressure growth 
persist iii 

With some of the materials questions coming 
under control, the observation trends are raising 
some serious questions. The decrease in T c follows 
dT c /dx ~ 1 K/% C. Ill wit h several reports of T c 
ks 21-23 K at s=0.2. |8f9ll0j However, with carbon 
adding one electron and the a bands of MgB2 hold- 
ing only ~0.11 holes per unit cell^ one might guess 
that the a bands would be filled (or nearly so, as 



the 7r bands can also accept carriers) and supercon- 
ductivity would have vanished. Not only is T c still 
robust at x=0.2, analysis of tunneling spectra in- 
dicates two-band superconductivity is retained; 14 ' 15 
whereas it might also seem likely that disorder scat- 
tering should have averaged out the gap to a single 
value. The upper critical field H C 2 initially increases 
strongly with C content as T c is depressed, reaching 
a maximumi*i2ii& of ~33-35 T. This critical field is 
substantially less than that for more disordered films 
(see above). 

Recent studies seem to agree that C can be intro- 
duced substitutional^ for B up to 10% replacement 
a;=0.10&i£ Tunneling spectroscopies indicate that 
Mg(Bi_ a; C K )2 remains a two gap superconductor^*^ 
(i.e. the gap anisotropy is not washed out by scatter- 
ing) for 10±2% C substitution for B. This question 
of how the anisotropy gets washed out has attracted 
much interest. Similarly, the mechanisms underly- 
ing the high critical fields are not understood. These 
questions are not the focus of the present study, al- 
though the knowledge gained from first principles 
calculations will be useful in the resolution of these 
topics. 

A few first principles calculations have been done 
for C substitution of B in MgB2. Pseudopo- 
tential calculations including relaxation in a 27 
unit cell supercelli^ confirmed that substitutional 
C is energetically favorable to interstitial C, and 
that the C-B bondlength is 5-6% shorter than the 
bulk B-B bondlength. A projector-augmented-wave 
calculationAS, reported a smaller relaxation, but the 
difference may be due to constraints related to su- 
percell size. A coherent potential approximation 
study of the disordered alloy using a Korringa-Kohn- 
Rostoker multiple scattering method in the atomic 
sphere approximation^ revealed relatively small ef- 
fects of disorder in the range 0< x < 0.3. The main 
effect of C substitution was reported to be the raising 
of the Fermi energy due to the additional carriers, 
but the filling of the crucial a hole band was not 
quantified. 

In this paper, we present a more detailed analysis 
of the changes in electronic structure, and the effec- 
tive a band doping, using both periodic supercells 
and disordered alloy calculations. Our study lays the 
groundwork for another issue in the superconductiv- 
ity of MgB 2 alloys that has received no direct at- 
tention from experimentalists: the effect of decreas- 
ing holes on the strength and character of electron- 
phonon coupling. Theoretical studies2i*22i2^ have 
predicted the following remarkable changes in doped 
MgB2i as the number of holes decreases (hence 
the Fermi wavevector of the cylindrical a Fermi 
surfaces decreases), (1) the substantial downward 
renormalization of the E2 ff bond-stretching modes 
with Q< 2kp does not change, (2) the coupling 
strengths of these already strongly coupled modes 
increases, and (3) the total coupling strength A, and 
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hence T c , remains unchanged. This scenario as- 
sumes two-dimensionality of the a bands which is 
only approximately true, and neglects the change 
of electron-phonon matrix elements which should 
be reasonable. This argument also neglects the 
fact that the applicability of conventional electron- 
phonon (Midgal-Eliashberg) theory, already some- 
what suspecfe2i2^ in MgB2, definitely becomes in- 
applicable as the Fermi energy decreases and ap- 
proaches more closely the bond-stretching mode en- 
ergy (65 meV). Item (2) definitely portends unusual 
dynamics related to the B-B stretch modes. We hope 
the current paper stimulates experimental investiga- 
tion into these questions. 



II. CALCULATIONAL METHODS 

Two methods of assessing the effects of substi- 
tutional C in MgB2 have been used in the work 
reported here. Our ordered impurity calculations 
have been performed using the full-potential lin- 
earized augmented plane wave code Wien2k2&, ap- 
plying the generalized gradient approximation^! to 
the exchange-correlation potential. The basis set is 
reduced by using the APW+lo method^, retaining 
the accuracy of the LAPW method. KK max was 
set to 7.00 which is a high quality basis for the s- 
p electron systems. The Brillouin zone integration 
for the supercell was carried out using 432 k points 
in the irreducible part of the zone. We have used 
the experimental results of Avdeev et <dm± for the 
lattice parameters of the 8.33% doped supercell and 
the 10% doped virtual crystal. 

For the disordered (randomly substituted) al- 
loy calculations, we have used the full-potential 
nonorthogonal local-orbital minimum-basis scheme 
(FPLO) 29 applied in the coherent potential approx- 
imation (CPA)iSi In the calculations, 1152 (45 x 
45 x 10) irreducible k points, and valence orbitals 
2s2p3s3p3d for Mg and 2s2p3d for B and C, were 
used. The implementation of the CPA in FPLO re- 
lies on the Blackman-Esterling-Berk theory^i that 
includes random off-diagonal matrix elements in the 
local orbital representation. 



III. ORDERED Mg(Bi_ x C a; )2, x = 0.0833 

Smaller concentrations of ordered C substitution 
for B require larger supercells. Our choice of super- 
cell was guided by the desire for a small enough C 
concentration to be relevant to address experimen- 
tal data, yet not so small as to make it difficult to 
distinguish the effect of C addition or to make the 
calculations unreasonably tedious. Since interlayer 
hopping is not an issue we will address, the supercell 
will involve enlargement of only a single B layer. The 
shape of supercell should have as small as possible 




FIG. 1: (Color online) Top view of the B-C layer of 
the MgoBnC supercell. Large (blue) circles denote Mg 
atoms, small light (yellow) circles denote boron, and 
small dark (red) circles denote carbon. The supercell 
boundary is outlined. 



TABLE I: Lattice parameter values used in the virtual 
crystal and supercell calculations. 



Lattice 


undoped 


10% doping 


8.33% doping 


Parameter 


MgB 2 


MgBi. 8 Co.^ 


MgeBnC 




(A) 


(A) 


(A) 


a 


3.083 


3.053 


3x3.053=9.159 


b 


3.083 


3.053 


2x3.053=6.106 


c 


3.521 


3.525 


3.525 



aspect ratio in order to maximize the separation of 
C atoms. The compromise we chose was the 8.33% 
doped supercell (Mg 6 B n C -> Mg(B 1 _ :c C :c )2, x = 
1/12 = 0.0833). Figure |T] shows the top view of the 
BuC layer in the 2x3 supercell. Each carbon atom 
has three first, five second, and three third nearest 
neighbor boron atoms. 



A. Band structure 

There are two effects of carbon substitution: a 
change in the average potential in the B-C layer, and 
the breaking of symmetry by C replacement of B in 
the supercell. In the ordered supercell the potential 
difference is kept explicit and symmetry breaking 
will be clear. In Sec. VI CPA will be used to probe 
these effects in a different manner. Another effect, 
B relaxation around the C impurities, has been ad- 
dressed previously to some extent, but it will not 
be studied here since the effect cannot be included 
within the CPA. 



FIG. 2: Band structure of MggBnC in the region of 
the Fermi level. The twofold degenerate a band in the 
undoped system along the F-A is split by 60 meV when 
1/12 of B is replaced by C. 



The band structure within 4 eV of the Fermi 
level (Ej?) is shown in Fig [21 The primary band 
in the understanding of MgB2 superconductivity is 
the (twofold degenerate) a band along the T-A di- 
rection. The results of supercell band-folding can be 
distinguished especially for the T-A direction, where 
(for example) the a band is replicated at T at -1.8 
eV, -2.3 eV, -2.6 eV, and even further lower energies. 

One of the most readily apparent effects of C dop- 
ing is the splitting of the a bands along T-A, by 
60 meV. It is helpful in understanding this splitting 
to obtain the atomic characters of each of the split 
bands, which will in any case involve only the B 
and C pa (p x ,p y ) states. The carbon p contribution 
is emphasized by enlarged symbols ("fatbands") in 
Fig. [3J The lower of the two bands has stronger C 
content (although this is not easy to distinguish in 
Fig. |3J). In addition, it has more first and second B 
neighbor character than does the upper band. The 
upper band contains more third neighbor B charac- 
ter, this being the B site farthest from the C atom 
and representative of the bulk material. The inter- 
pretation of this splitting is that lower band has 
been pulled down due to the stronger potential of 
the C atom compared to that of B. The 60 meV 
splitting for £=1/12 provides (assuming a linear ef- 
fect in this concentration range) an energy scale for 
CT-band broadening 7 ~ 7 meV/% C content. 




FIG. 3: Bands of MgeBnC near the Fermi level, with 
symbol size proportional to the C 2p character. The a 
bands, split by 60 meV (hardly visible in this figure) lie 
just above Ef at Y and disperse upward toward the A 
point. 



B. Density of States 

Figure0]shows the total and atom-projected den- 
sity of states (DOS) of the 8.33% doped system. The 
C and B DOS are similar at and above Ep, but the 
C DOS is somewhat lower in the interval within 2 
eV below E F . The C DOS is larger than that of 
B in the -8 eV to -3 eV region. Based on charge 
within spheres of 1.65 A radius we calculated a charge 
transfer of about 0.095e _ to the carbon from the 
three first nearest neighbor boron atoms, with other 
B sites showing negligible change in charge. This 
charge transfer also reflects the stronger potential 
and larger electronegativity of C, and provides sub- 
stitutional C with definite anionic character. 



IV. 10% DOPING 

We have also performed a virtual crystal 
calculation for the 10% doped 1=0.10 system 
Mg(B .9Co.i)2- The DOS and band structure are 
shown in Figs. [S] and |H| By aligning the highest 
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FIG. 4: Density of states, total (bottom) and decom- 
posed into B and C contributions on a per-atom basis, 
for the ordered x=l/12 model. The primary B - C dif- 
ference is that the C DOS is lower in the -2 eV to E_f 
region, but higher in the lower region -8 eV to -2 eV. 
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FIG. 5: Density of States of Mg(Bo.gCo.i) 2 . The Fermi 
level moves up by 0.3eV due to doping 



peaks in the DOS at -2 eV and +6 eV (they can be 
aligned simultaneously) , we establish that the differ- 
ences are (1) an increase in the occupied bandwidth 
(at the MgB 2 band filling) from 12.4 eV to 12.5 eV 
(~1%), and (2) the raising of the Fermi level by 0.3 
eV to accommodate the extra electrons. Otherwise 
there is very little difference in the two densities of 
states. This raising of Ep will be compared to the 
supercell result in the next section. At this band fill- 
ing, the VCA gives the a band edge at T precisely 
at Fjp as shown in Fig. HJ1 At the point A it is still 
about 0.3 eV above E F , meaning that both cylindri- 
cal Fermi surfaces have shrunk to a point at their 
"waists." 
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FIG. 6: Band structure of Mg(Bi_ a; C :c )2, x = 0.10, in 
the virtual crystal approximation. The a bands at Y he 
precisely at Ef, so the Fermi cylinder radii at F vanish, 
corresponding to a topological transition. 



V. a HOLE CONCENTRATION 

It is accepted that superconductivity in MgB 2 
arises from hole-doping of the a bonding hole 
states due to the intrinsic chemistrjii^ of MgB2. 
The observed decrease in electron-phonon coupling 
strength, and hence T c upon C addition, makes the 
change in hole concentration one of the major points 
of interest. We have calculated the hole concentra- 
tion of both the 8.33% supercell system and the 
10% virtual system from the volume enclosed by 
the Fermi surface. Since the radii of the cylindri- 
cal Fermi surfaces are only slightly different along 
the r — > M and the r — > K directions, we calcu- 
lated the average basal area by considering the value 
along the r — » M direction only. Also, we assumed 
a sinusoidal dependence of the Fermi surface along 
the c axis, allowing analytic evaluation of the Fermi 
surface volumes i.e. the hole concentrations. The 
results are presented in Table II. In the process of 
calculating the number of holes, electronic structure 
reveals that both the a band Fermi surfaces are still 
intact, consistent with the two-band superconduc- 
tivity with substantial T c as seen in experiments^ 
Also, we interpolated the number of holes for the 
10% doped system (virtual crystal) to calculate the 
value for the 8.33% doped system. Supercell calcu- 
lation gave 0.070 holes/cell while the extrapolation 
gave only 0.057 holes/cell. This indicates that the 
virtual crystal approximation is not very reliable for 
substitutional carbon, and that C cannot be thought 



TABLE II: Calculated number of a holes for the various 
C concentrations discussed in the text. 





JNo. 01 holes 




(holes/cell) 


MgB 2 


0.11 


8.33% doping 




Mg6BnC supercell 


0.070 


10% doping 




MgBi.gCo.2 virtual crystal 


0.0463 


8.33% doping 




MgB1.g33C0.i67 extrapolation 


0.057 



of "boron + an electron." 



VI. COHERENT POTENTIAL 
APPROXIMATION RESULTS 

For further comparison and to assess the effects 
of disorder, we have performed CPA studies of this 
system. We have carried out CPA calculations (1) 
at .x=0.0001 to provide a "perfect crystal" reference 
for evaluation of the various algorithms in the CPA 
code, (2) at x = 0.0833 for most direct comparison 
to the x=l/12 (ordered) supercell calculations de- 
scribed in the previous sections, and (3) at x= 0.10 
and 0.20 as representative of the system toward the 
upper range of achievable C substitution. The spec- 
tral function over the whole energy range of interest, 
plotted as a "smeared" band structure, is shown for 
the 1=0.20 case in Fig. [7| The primary points of 
interest are the filling of the a band hole states, and 
the broadening (and potentially splitting) of bands. 

The full energy region is presented in Fig. {7\ for 
£=0.20 (where broadening is more easily seen) to il- 
lustrate the strong wavevector (k) and energy (E) 
dependence of the broadening of the spectral func- 
tion. The largest disorder occurs in the 2s region in 
the lower valence band, where the C 2s state is no- 
ticeably lower in energy leading to increased smear- 
ing. The other region of large disorder is from +3 
eV upward, but some states remain comparatively 
sharp. The flat band at 5 eV along T-M and A- 
L, with strong n character, is the lowest conduc- 
tion band that is strongly affected by the disorder. 
The bands around the Fermi level, whether a or tt, 
are among those less affected by the chemical dis- 
order. As expected from the virtual crystal results 
for ir=0.10 from the previous section, the a band 
holes are completely filled at £=0.20, as revealed by 
the flat band along the T-A line lying entirely below 
Fip- The band filling behavior versus x is quantified 
below. 

In the valence (occupied) bands, disorder broad- 
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FIG. 7: CPA spectral density in the full valence- 
conduction band region, for x=0.20, plotted as a broad- 
ened band structure. Disorder broadening is largest be- 
low -8 eV in the valence bands and above 4 eV in the 
conduction bands. 



cning is large in the lower s-band region -8 eV to 
-14 eV below Ep, and somewhat less so where the 
(primarily a) band in the -3 to -4 eV range gets flat 
around the zone edge M and L points. Similar broad- 
ening does not occur in the same band at the (more 
distant from T) zone edge K and H points, where 
the band lies at -6 eV and is much less flat. The top 
of the a bonding bands along Y- A are comparatively 
sharp; the width is quantified below. 

In the conduction bands disorder broadening be- 
comes more prevalent. The antibonding a* bands 
(flat along T-A at +6 eV) are much broader than 
their bonding counterparts (below Ep), and this 
band along T-M and A-L with Mg 3s character be- 
comes exceedingly diffuse. The 7r bands show strong 
k-dependence of the broadening, beginning at 2-3 eV 
around the M point and becoming wider in the 6-10 
eV range (and above, not shown in the figure). 



Carbon concentration dependence. The x depen- 
dence of the broadening can be seen in Fig. [S] for 
.x=0.0001 (the CPA equivalent of x=0 MgB 2 ), for 
x=0.10, and for £=0.20 (the latter may not be ex- 
perimentally accessible). Since in each case Ep is 
set to zero, band filling appears as downward shifts 
of the bands, by roughly 5Ep — -0.4 eV for x=0.10 
and SEp = -1.0 eV for x=0.20. For x=0.0833 (not 
shown) the a bands lie at the same position (with re- 
spect to Ep) as for the ordered x=l/12 case shown 
in Fig. [5] Hence the degree of band filling in the 
CPA results is the same as for the ordered £=1/12 
case that was analyzed in the previous section. Al- 
though we will not dwell on it, it can be noticed that 
the band shift is not entirely rigid: the two valence 
bands at the L point split apart as well as broaden 
with increasing C concentration. 

The other important aspect of the CPA bands of 
Fig. |2I is the broadening that increases with carbon 




FIG. 8: CPA spectral density for £=0.0001 (top panel), 
i=0.10 (middle) and £=0.20 (bottom panel). The main 
features are (i) the "rising" Fermi level (with respect to 
the a band along T-A, say), and (ii) the increase in the 
disorder broadening of the bands. Since the numerical 
algorithms cannot reproduce the <5-function bands for 
x — > 0, the £=0.0001 case is included as a reference for 
the algorithmic contribution to the width. 



concentration. Within the resolution of widths that 
we are able to extract, the widths of the a bands 
can be taken to be proportional to x, but differ con- 
siderably between T and A (being about 50% wider 
at r) . The full width at half maximum of the spec- 
tral density of the a bands, averaged between T and 
A, is about 7 «0.21 eV (to perhaps 10% accuracy). 
This width corresponds to a width in wavevector 
given by 7 = vpSk, and therefore a mean free path 
of lp = 2irvp/j. The 60 meV splitting of the a 
bands (see Sec. III. A) for the x = 1/12 ordered su- 
percell implies that the 0.21 eV mean width can be 



interpreted as 0.05 eV from the B/C on-site energy 
difference, and thus 0.16 eV from disorder itself. 

Due to the anisotropy of the Fermi velocity vp, 
the mean free path may vary considerably over the 
Fermi surface at 1=0.10. At T the cylinder radius 
has shrunk to a point, and the very small z com- 
ponent of vp (vanishing at T) suggests a very small 
mean free path in the z direction of the order of the 
layer spacing c for x=0.10. 



VII. DISCUSSION AND SUMMARY 

In the band structure of the 8.33% doped super- 
cell system, the hole a-band Fermi surfaces along 
the r — > A direction are still present, and this de- 
gree of band filling is reproduced by the CPA. Even 
for a;=0.10 the virtual crystal picture (which overes- 
timates the rate of band filling by C) leaves a Fermi 
surfaces, just beginning to be pinched off. The CPA 
calculations (Fig. [SJ show the a band holes begin 
to disappear rapidly for x > 0.10. Qualitatively this 
filling is consistent with most experimental reports. 

The questions, and experimental probes, should 
now be focused on the fact that, as the a band 
fills, there will be very strong deviation from 'busi- 
ness as usual' in the coupled electron-phonon sys- 
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tern. The strength of coupling of bond-stretching 
modes with Q< 2kp continues to increase, and con- 
ventional Migdal-Eliashberg theory ceases to apply. 
The dynamics of these ultra-strongly coupled modes 
is unexplored, with their peculiar character being 
signaled by the divergence of their linewidth (at least 
within Migdal-Eliashberg theory). The limiting be- 
havior should not revert to the widely studied po- 
laron limit, however, as there remain the background 
7r electrons, which are weakly coupled to vibrations 
but provide full metallic conductivity and screen- 
ing to the system. The present study reveals that 
C substitution for B provides a similarly favorable 
system to Al substitution for Mg for studying this 
evolution. Recent reports indicate that Sc substitu- 
tion (Mgi_ :c Sc a; B2) may also provide^ another such 
system for study. 
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